Analysis of the phase transition for the Ising model on the frustrated square lattice 
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We analyze the phase transition of the frustrated J1-J2 Ising model with antiferromagnetic nearest- 
and strong next-nearest neighbor interactions on the square lattice. Using extensive Monte Carlo 
simulations we show that the nature of the phase transition for 1/2 < J2/J1 < 1 is not of the weakly 
universal type - as commonly believed - but we conclude from the clearly doubly peaked structure 
of the energy histograms that the transition is of weak first order. Motivated by these results, we 
analyze the phase transitions via field-theoretic methods; i.e., we calculate the central charge of the 
underlying field theory via transfer-matrix techniques and present, furthermore, a field-theoretic 
discussion on the phase-transition behavior of the model. Starting from the conformally invariant 
fixed point of two decoupled critical Ising models (Ji = 0), we calculate the effect of the nearest 
neighbor coupling term perturbatively using operator product expansions. As an efi'ective action we 
obtain the Ashkin- Teller model. 

PACS numbers: 64.60.De, 75.10.Hk, 05.70.Jk, 75.40.Mg 



I. INTRODUCTION 

The simple Ising model which adds up the interac- 
tions of two-state variables on a D-dimensional lattice 
has served as a pioneer in many physical problems, es- 
pecially in statistical mechanics and solid-state physics.^ 
It was one of the first models to mimic the magnetic 
exchange interactions in condensed-matter theory, and 
the behavior of phase transitions was studied extensively 
for this model in different dimensions and by numerous 
approaches. Nevertheless there arc still open questions 
concerning this model, in particular in two and three di- 
mensions. In two dimensions the model with nearest- 
neighbor (NN) interactions undergoes an ordering pro- 
cess at a finite temperature which is well understood and 
establishes the Ising universality class for second-order 
phase transitions.^ However, for frustrating interactions 
the phase diagram of the model becomes richer and the 
physics of the occurring phase transitions becomes more 
complicated. 

In the present work we focus in particular on the 
phase transition from the high-temperature paramag- 
netic phase into an antiferromagnetic coUincar phase 
that is favored by strong additional interactions on next- 
nearest-neighbor (NNN) bonds, i.e., J2 > Ji/2. This 
transition also attracted a lot of interest in the past. 
In the late 1970s first rcnormalization-group calcula- 
tions and Monte Carlo (MC) simulations for the Ising 
model with frustrating interactions were performed by 
Nightingale'^ and Swendsen and Krinsky^ and later on 
by Oitmaa^ and Landau and Binder. They assumed a 
continuous phase transition and computed for this par- 
ticular model transition temperatures and critical expo- 
nents. Throughout the 1980s it was commonly accepted 
that the exponents are weakly universal but vary for dif- 
ferent degrees of frustration; i.e., only the scaling re- 



lations for the exponents are fulfilled but the absolute 
values are not universal. However, a continuous phase 
transition with non-universal exponents is only possible 
if the central charge of the underlying conformal field the- 
ory (CFT) is c > 1}^ Meanwhile theories with discrete 
c < 1 define universality classes with universal exponents 
such as the two-dimensional Ising model. ^^'^^ Since the 
present system can be described by two copies of Ising 
models in one part of the phase diagram, the universal- 
ity of the phase transition was under debate. In 1993 
Lopez et al. presented a mean-field calculation for the 
model^'^'^^ where they find a first-order transition for a 
finite parameter region of 0.5 < J2/ Ji ^ l-l- Recent MC 
simulations by the group of Malakis et al. contradict 
this scenario at least for the value J2I J\ = 1-^^ On the 
other hand, our MC results in Refs. 16 and 17 strengthen 
the scenario of a first-order transition for small values of 
0.5 < J2/J1 < 0.7. 

Here, we will show that the first-order scenario is valid 
up to J2 < 0.9 J\. However the MC simulations do not 
give a conclusive picture for larger values of Ji due to in- 
creasing length scales. Hence, it was necessary to apply 
further techniques to understand the nature of the phase 
transition for these parameters. To compute the central 
charge of the model we performed a finite-size analysis 
of the free energy which was calculated via transfer ma- 
trix techniques, but we can only get reasonable results 
for large J2/J1 > 1. In a last step we start from the 
limit of two decoupled Ising models (see Fig. 1 below) 
with antiferromagnetic coupling J2 and add perturba- 
tively an antiferromagnetic nearest-neighbor interaction 
Ji between the two copies. In second order wc arrive 
at an Ashkin- Teller model which is in agreement with a 
scenario of non-universal exponents. 

The paper is structured as follows: we will present the 
model in detail in the following section before review- 
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FIG. 1. (Color online) CoUinear phase of the Ji — J2 square 
lattice. Red dots stand for up spins, and yellow dots stand 
for down spins. The two copies A and B of the Ising model 
with magnetic couplings J2 and lattice spacing a are, respec- 
tively, represented with dashed clear and dotted dark blue 
lines, while the black thin lines correspond to the Ji square 
lattice. The shaded area represents the unit cell used to de- 
rive the continuum limit (Sec. IV B). The coordinates are 
indicated with respect to the x- and y axis of the A sublat- 
tice. 



ing the MC simulations in Ref. 16 and new results that 
we will present in Sec. III. We will also show that the 
length scales, which are needed to see the first-order fea- 
tures, are growing with J2 and are not accessible any 
more with MC simulations for J2 > Ji. The same scal- 
ing problems occur for the transfer matrix calculations in 
Sec. IV A, which are used to determine the central charge 
of the underlying field theory in the limit of J2 \ Ji- In 
Sec. IV we derive the conformal field theory for the case 
of two independent Ising models which interact via a per- 
turbation caused by the nearest-neighbor coupling Ji. 



II. MODEL 



The lattice model is described by the Hamiltonian 
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where the sums run over all nearest-neighbor and next- 
nearest-neighbor interactions on a N = L x L square 
lattice (see Fig. 1). The energy on each bond is given by 
the product of the adjacent classical Ising spins Si = ±1 
and the corresponding Ji which arc both chosen antifer- 
romagnctic. Thus, a configuration which yields minimal 
energy for all bonds does not exist for finite Ji and, hence. 



the model is frustrated. For small competing interac- 
tions J2 < Ji/2 on the diagonals of the square lattice, 
the model undergoes a phase transition from the para- 
magnetic phase into a Neel ordered configuration at a 
critical temperature Tc{J2) which depends on the frus- 
tration (compare Fig. 2). This phase transition is contin- 
uous and the scaling exponents arc the same as for the 
unfrustratcd square lattice Ising model. ^^^^ If Ji = 2 J2 
the critical temperature is suppressed to zero and the 
ground state is degenerate of order 2^ }^'^'^ For smaller 
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FIG. 2. (Color online) Critical temperatures for the phase 
transition from the paramagnetic into the magnetically or- 
dered phase over the strength of frustration J1/J2. We 
adapted the energy scale of the temperature for J1/J2 < 2 
to J2 and for J1/J2 > 2 to Ji. Note that the frustration is 
also given in units of J2/J1 on the upper x-axis. 

Ji the ground state is a coUincar antiferromagnet where 
lines of parallel spins are coupled anti-parallel. 



III. MONTE CARLO SIMULATION 

Based on the results we have presented in Rcfs. 16 
and 17 we performed further MC simulations on larger 
lattices and for new parameters. We used a single-spin 
Metropolis MC update^* with an additional line update. 
For an optimized thermalization process and high-data 
quality we also implemented an exchange Monte Carlo 
update^^"^^ and ran the simulations on large-scale clus- 
ters using OpenMP and MPI. The statistical errors of the 
data were obtained by multiple independent MC runs. 

The phase diagram [Tc over J1/J2) showing both old 
data and new transition temperatures up to J2 = 10 Ji 
is given in Fig. 2. Additionally we looked at the char- 
acter of the finite-temperature phase transitions and cal- 
culated critical exponents via finite-size scaling from our 
MC data. Comparing with old results from Landau and 
Bindcr^^^ we found a discrepancy between their values 
and ours for J2 > 0.5 Ji. To have a closer look at the 
nature of the phase transition for this part of the phase 
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FIG. 3. (Color online) Energy histograms (solid black lines) 
for two values of J2 = 0.8 Ji (left) and J2 = 0.9 Ji (right) 
for a L = 1000 (left) and L = 2000 (right) lattice and the 
reweighted histograms for slightly lower temperatures (dashed 
red lines). The two-peaked structure emphasizes the first- 
order character of the phase transition. 



diagram we recorded energy histograms at discrete tem- 
perature steps. These are plotted for J2 < 0.7 Ji in Fig. 4 
of Ref. 17 and show a clear two-peaked structure and 
therefore prove the first-order character of the phase tran- 
sition. Since Lopez et at claim the first-order transition 
scenario to be valid up to J2 ~ 1.1 Ji,^^'^^ we recorded, 
for the present work, histograms also for larger values of 
J2 . In Fig. 3 we show the recorded histograms of the MC 
simulations (as solid black lines) for a 2000 x 2000 lattice 
at J2 = 0.9 Ji and a 1000 x 1000 lattice at J2 = 0.8 Ji. 
Note that for a small step in the parameter space (0.1 Ji) 
it is already necessary to double the linear size of the 
simulated lattice to achieve a similar resolution for the 
recorded histogram. The shapes of these histograms show 
a strong deviation from the almost Gaussian shape that 
is expected for a continuous phase transition, yet for the 
given temperatures the structure is not symmetric. For 
this reason, we also present reweighted histograms (as 
dashed red lines) for slightly lower temperatures, i.e., the 
size-dependent transition temperatures. The histograms 
are thus shifted to lower energies and exhibiting a more 
distinct and symmetric two-peak structure. The simula- 
tion and recording of a new histogram for this tempera- 
ture would have been too time consuming and, hence, we 
applied the standard reweighting technique^^'^'^ to prove 
the first-order character of the phase transition. 

For larger J2 the resolution of a double peaked profile 
in the histograms is not feasible, probably due to the 
growing crossover scales. 

In addition we present measurements of the spin-spin 
correlation function {StSj). In Fig. 4 we show separately 
the correlation functions {SiSj) for i and j being part of 
the same sublattice (top) and being in different sublat- 
tices (middle) for the value Ji = J2 . For the correlations 
inside the same sublattice - in addition shown in a dou- 
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FIG. 4. (Color online) Calculations for the spin-spin corre- 
lation functions (SiSj) at Ji = J2 for a 100 x 100 periodic 
lattice along one direction. Thus maximal distance is given 
by i — J = 50. Top: Correlations inside the same sublattice - 
i — j even - for different temperatures around Tc- Middle: 
Correlations between spins on different sublattice sites - i — j 
odd - in a larger scale. For all temperatures the correlations 
decay rapidly and go to zero. Hence, no long-range corre- 
lation is observable between the two sublattices. Botttom: 
Doubly logarithmic plot of the correlation functions inside the 
same sublattice for three exemplary temperatures and their 
related fits. 



bly logarithmic scale (Fig. 4 bottom) ~ we observe the 
behavior of a second-order phase transition, i.e., expo- 
nential decay for T > Tq, algebraic decay for T « Tc and 
fast saturation toward a constant finite value for T < Tc- 
We also extracted the critical exponent 77 which describes 
the scaling of the correlation in the vicinity of the crit- 
ical temperature and obtained ?/ = 0.20(1). This is far 
away from the Ising value of 77/ = 0.25. But for smaller 
values of Ji / J2 the critical behavior becomes more Ising- 
like and we obtain values of 77 = 0.25(1) for J1/J2 < 0.2 
(not shown). In the middle part of Fig. 4 we show the 
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correlations between spins on different sublattices which 
decay for all temperatures quicker than correlations in 
the same sublattice - note the different scales for the two 
upper panels of Fig. 4. For T < Tc the correlations drop 
to zero. We observed this behavior for all J1/J2 < 2, i.e., 
in the region where the system undergoes a phase tran- 
sition into the coUinear phase. But the decay in (SiSj) 
becomes slower for increasing J\l Ji. 

In conclusion the MC data yield a clear picture only 
for 0.5 Ji < J2 < J\ where a first-order phase transition 
scenario is established by the doubly peaked structure of 
the energy histograms. For larger values of J 2 the ana- 
lysis of the correlation functions indicates a decoupling 
of the two sublattices and a continuous phase transition. 
However, a detailed examination of the scaling behavior 
at the critical temperature and a reliable calculation of 
critical exponents are hampered by large crossover scales. 

IV. CONFORMAL FIELD THEORY 

In the limit Ji = 0, the system is exactly described 
by two decoupled two-dimensional Ising models on sub- 
lattices A and B (see Fig. 1), and the critical behavior 
can therefore be tackled by conformal field theory. ^^'^"'■^^ 
A CFT is characterized by a constant c called the cen- 
tral charge (or conformal anomaly). Various CFTs with 
c < 1 were identified as statistical models at their crit- 
ical point. In particular, the CFT with c = 1/2, which 
was known to correspond to the massless free Majorana 
fermion, was identified as the critical two-dimensional 
Ising model. ^^'^^ Physically, the central charge character- 
izes the short-distance behavior of the theory and can be 
seen as a measure of the number of degrees of freedom of 
the system. Central charges of decoupled systems add up 
and, therefore, in the limit Ji = where the two Ising 
models A and B are independent, it is c = 1. At this 
point we want to emphasize that every unitary theory 
with c < 1 yields a universality class with constant uni- 
versal exponents. Hence, to gain a phase transition with 
varying critical exponents is only possible for a CFT with 

C> 1.26 

A. Transfer Matrix 

In the 1980s it was shown that, for a lattice model, the 
free energy / per site of a cylinder of infinite length and 
finite circumference L shows a finite-size scaling in 1/L 
with a proportionality factor depending on the central 
charge of the corresponding field theory.^^'^* 

/ = a - c"^^ + OiL-^) , a = const. (2) 

To calculate the free energy of our model, MC simula- 
tions are not suitable since the entropy is not available 
for them. Therefore, we chose to implement a transfer- 
matrix algorithm. 2'^^ We were able to calculate the free 
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TABLE I. Central charge of the underlying field theory for 
dilTerent J2/J1 calculated with a transfer matrix computa- 
tion of the free energy. The width of the computed systems 
satisfies L < 22. We also include in parentheses c-values for 
Ji/2 < J2 < Ji where the phase transition is of weak first 
order according to the MC analysis. 

energy for systems of size Lx B, where the circumference 
of the cylinder is limited to L < 22 because of exponen- 
tial growth of computational effort and the length can be 
chosen easily up to i? = 10000 (linear scale). A finite-size 
scaling of the free energy yields an estimate of the central 
charge. We present our results for different J2/J1 in Ta- 
ble I. For large values of J2 the central charge seems to 
converge to the value of two independent Ising models, 
c = 1. However, for decreasing J2 they do not converge 
and we find c > 1, which is not allowed by the c theo- 
rem which states that the central charge can not increase 
under the influence of a renormalization-group transfor- 
mation if the corresponding field theory is critical. In- 
deed the MC data indicates a weak first-order transition 
for J2 = 0.6 Ji^s and J2 = 0.8 Ji (Fig. 3 left-hand side) 
accompanied by large crossover scales. Thus, it is not 
surprising that this weak first-order transition is not de- 
tected by the transfer matrix computations for cylinders 
with circumferences i < 22 such that the corresponding 
results for c are not meaningful. 



B. Continuum limit 

In order to get more information on the nature of the 
phase transition for Ji ^ we now derive our model 
in the continuum limit. We start from the conformally 
invariant fixed point Ji = (see Fig. 1) and then add 
perturbatively a coupling oc Ji between the two decou- 
pled Ising models A and B. 

In a first step, the spin variables on discrete lattice 
sites are replaced by continuous fields which resemble the 
underlying Neel order of the two decoupled antiferromag- 
netic Ising models. Thus, the staggered spin variables Sj 
{I = A,B sublattice) have to be transformed into smooth 
variables defined by 

CT/(m, n) oc (-l)™+"S/(m, n) , (3) 

where (m, n) are the lattice coordinates. Note that the 
transformation Eq. (3) is based on a specific choice of 
gauge but this does not affect the macroscopic proper- 
ties of the system. The choice of the unit cell is shown in 
Fig. 1 and the coordinate system is rotated by an angle 
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of 7r/4; i.e., the axes point along the next- nearest neigh- 
bor bonds. In a next step the sum of the Hamihonian 
equation (1) is converted into a two-dimensional integral 
where the values of <7i{x,y) at the limits (— cx), cx)) are 
equal due to the periodic boundary conditions imposed 
on our model: 



E 



a: lattice spacing (4) 



A Taylor expansion up to second order on the aj = 
ai{x,y) fields is given by: 

a I {x + m a, y + n a) = aj + a [mSx + ndy] cr/ 

+ y [m'd^, + n'dl + 2mndl^] + 0{d^a,) . (5) 

By summing up all nearest- neighbor interactions ( Ji) ap- 
pearing in the chosen unit cell, the products cr a<^b and 
aidx,y(J.j {I 7^ J) are canceled out due to frustration.-^^ 
This feature distinguishes the present model from the 
two-layer Ising model, where the highly relevant ctacb 
coupling survives. '^^''^'^ The Ji interaction is then given 
by: 



Hint — ^A*' 



'Ji jj dxdy{dx(JAdy<7B + dyUAdx^B^ 



Note that partial derivatives are understood to act only 
on the subsequent operator and is a constant factor. 

For further calculations, it is convenient to rewrite the 
interaction [Eq. (6)] in complex coordinates z = x + iy 
and z = .T — iy: 

nint^-lf^^Jl JJdzdz(^0i{z,z)~02{z,z)^ (7) 

Oi{z, z) = dz<7Adz<JB , 02{z, z) = dgaAdgas ■ (8) 

The CT/ fields of the theory have conformal dimensions 
{ha,K) = (1/16, 1/16). "'24-26 One defines the scaling 
dimension A ^ h + h and the conformal spin s = h — h. 
The full model is then described by the action 

A^A\+A%+T jj dzdz [eA{z, z) +eB{z, z)^ (9) 

+ g jjdzdz(0i{z,z)-02{z,z)y g(xJi, (10) 

where are the fixed-point actions of the Ising models 
A and B. Furthermore £/ are the usual thermal operators 
of the conformal field theory on the two-dimensional Ising 
model with {he, he) = (1/2,1/2) and the corresponding 
mass T oc (T - T^'="). 

The operators Oi and O2 arc kept, despite being highly 
irrelevant with a scaling dimension A = 9/4, due to 
their non-zero conformal spins s = ±2. The presence of 
such chiral terms was previously reported in other frus- 
trated systems such as the anisotropic square lattice, 
the checkerboard^^ and the Kagome'^^ lattices. These 
twist terms^"^ are known to be likely to generate relevant 
or marginally relevant terms at higher order. '^^''^^ 



Before calculating higher orders of the perturbative in- 
teraction, we want to briefly discuss the underlying sym- 
metries of the model and the consequences for the contin- 
uous field theory. The Hamiltonian Eq. (1) is invariant 
under translations by multiples of the lattice spacing a, 
which is preserved for the integral form of the field the- 
ory because of the same periodic boundary conditions 
imposed on the integrals. Furthermore, the model is 
symmetric under rotations by an angle of tt/2 and in- 
versions along the axes {x,y = 0), {x = 0,y), {x,y = x) 
and {x, y ~ —x) (corresponding to the diagonals and the 
vertical and horizontal line through the origin in Fig. 1 
residing on an A site of the lattice). Note that the ro- 
tations and inversions are not independent. In addition 
the model is symmetric under the exchange of the two 
sublattices, which is equivalent to the translation by one 
lattice spacing of the original lattice. The inversion, e.g., 
of the axis {x,y = 0) — >■ {—x,y = 0) is given in the 
complex coordinates by setting: 



z — >■ —z, z — > 



dz -> -dz, dz -5, 



tT^(z, z) aA{—z, —z) and 
0-5(2, z) -asi-z, -z) 



(11) 



Since the cr/ fields are completely symmetric in z and 
z and we integrate over the whole complex plane for 
both variables, the crucial point of this transformation is 
the change of the partial derivatives and the additional 
sign that occurs for the B sublattice. Thus, operators 
that contribute to the continuous field theory have to 
be either quadratic in the sublattice fields or contain an 
asymmetric contribution of derivatives and fields resid- 
ing on different sublattices, as those in Eq. (8). Further- 
more, the inversion (x, y — —x) — > {—x, y = x) rotates 
the partial derivatives onto the imaginary axis, which 
ensures that by symmetry only operators containing an 
even number of partial derivatives are allowed. Following 
these symmetry arguments, we can discuss the appear- 
ance of certain operators in higher orders. A highly rele- 
vant (A = 1/4) spin-spin coupling aA{z,z)aB{z,z) such 
as the one appearing in the two-layer Ising model^^.sa jg 
not allowed, whereas combinations of energy operators 
ej are allowed since they transform like products of spin 
fields on the same sublattice. 



C. Operator product expansion 

Higher-order perturbations are calculated via the stan- 
dard operator product expansion (OPE). ^^'^4, 25 'pj-^jg Qp_ 
eration allows us to replace, inside a correlation function, 
a product of two operators by a combination of scaling 
operators allowed by the theory. This is meant to close 
the renormalization-group equations in the operator al- 
gebra of the model before discarding irrelevant pertur- 
bations. For the two-dimensional Ising model, the field 
content of the product of two fields is encoded in the 
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fusion rules: 



M ['Jj]=Sjj{[l] + [ej]) 
N [ej]=Sij[l] 



(12) 



where 1 is the identity operator present in ah CFTs and 
I, J represent the two sublatticcs A and B. Using the 
general normalized form of an OPE and the fusion rules, 
one obtains:^^ 



ai{z,z)aj{w,w) 



\z—w 



'''■rT7i+'-¥\^-M'^^ei{w,w) 



+ - - w)3/8a„ej(u;, w) + H.c. . (13) 

Note that we kept marginal terms with non-zero confor- 
mal spin to be consistent with previous comments about 

I 



the relation of the eventual importance of chiral terms in 
frustrated systems. 

By calculating the OPE between the ©1^2 operators 
of Eq. (8) one generates higher-order terms. At second 
order this generates terms [e = £{w, w)]: 



CX £/, OC 9{„,,iB}£/, OC £/£j (14) 

« £/i9{u,_iB}£j and oc (9{u,,ts}£/i9{u;,ts}e.7 (li^J) ■ 



(15) 



Following the calculations presented in Appendix we 



get: 
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[Oi{z,z) - 02iz,z)][0i{w,w) ~ 02(10, w)] = ao (z ~ w) 4 (z - w) 4 + H.c. 
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4096' 



ai 



a2 



(z — w) 1- [z — w) + H.c. 



(z — w) [z — w) + H.c. 



4096 

81 
8192 



\z-w\ 2 (. (eA + SB) 
z~w\ '^ > (£a£b) 



— (^a3{z — w) i{z ~ w) -1 + 0:4(2 — ui) i {z — w) -i + a-i{z — w) ^ [z — w) 1 ^ {dw£A + dw£B) + H.c. 

/ 9 3 1 5 7 13 \ 

+ f 05(2 — w)^ 3(2 — u)) 4 — af,{z — w)^ ^ [z - w)^ ^ + 05(2 — w) 4 (z - \ {dw£A£B + £a<9„,£b) + H.c. 

+ ^ /3fc(z - {Z - wf- {0{dl^^^^£A£B)) , Vk^tk 



(16) 

(17) 

(18) 
(19) 
(20) 



The coefficients ai and j3k are some rational constants, 
and the terms in Eq. (19) are irrelevant but could - as 
shown before - produce again relevant terms in higher or- 
der. All terms in Eq. (20) are also irrelevant and contain 
second derivatives that will produce only highly irrele- 
vant terms (A = 4) in higher orders. 

The third-order oc terms are obtained by multiply- 
ing the above operators with the original perturbation 
Eq. (8). One needs the following OPE: 

Cr/(z, z)£j(w, w) — -^1^ — w\^^(Jl{w, w) 

-h%^(z -w)^{z- w)-^d^ai{w, w) + H.c. (21) 

It yields only three different types of operators that are 
primary and secondary operators from the spin family: 

CX aiaj , CX (Tid{y,^y:,}aj (22) 
and CX d{.^^^,}aid{.^^^}aj . (23) 

Calculating the third-order perturbation by multiplying 
every operator from the second order with Oi — O2 yields 
the prefactors for these spin operators. However, for each 
of the relevant operators in Eq. (22) the sum of all pref- 
actors gives exactly zero, which is in agreement with the 



symmetry considerations presented in Sec. IV B. Thus, 
the third order does not give any new operators since the 
irrelevant terms in Eq. (23) are the same as in the first 
order of our perturbation. Thus, we have closed the oper- 
ator algebra and are left only with the terms of Eq. (14), 
namely, the thermal operators £a + £b and two marginal 
operators, £^£5 and (9{.u,,ts}£A + d{w^.u,}£B- 

The next step is the integration of the prefactors given 
in Eqs. (16), (17) and (18) which depend on the product 
of {z — wY and (z — w)*. Since for the second-order terms 
in Eq. (14) we have to integrate over all four variables, 
z, z, w, and w, with the only constraint being |z— luj > a, 
all prefactors with exponents v ^ t and \v — t\ — n (with 
71 G N) will be zero due to the phase integration from 
to 2tt. In detail one can carry out first the integral 
over z and z, which leaves the energy fields untouched; 
remembering that z and z are complex conjugates, one 
can rewrite (z — w)"" {z — wY = |z — w|''"*'*e"^^"~*-' , whereby 
4> ~ arg(z — w). 

This integration cancels the marginal chiral terms dw£i 
and d^£j, again in agreement with the symmetries of 
the model. Thus, the perturbation is given by the pure 
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energy terms and the energy-energy coupling 



4096 
8192 



dzdzdwdw \z — w\ '^{ea + sb) (24) 



\z—^a I >a 



/ dzdzdwdw \z — w\ '^Ea^b- (25) 

J \z—w I >a 



The signs of these operators stem from the squares of 
the coupHng given in Eq. (7) and the derived signs of 
the perturbative calculations given in Eqs. (16) and (17), 
respectively. The pure energy terms in Eq. (24) and the 
mass term in Eq. (9) are necessary to tune the model 
onto criticality. We find the new critical temperature up 
to second order in perturbation theory 



Tc{Ji) = Tc(0) 



(26) 



which appears to be in agreement with the behavior of 
the critical line on the left-hand side of the phase diagram 
presented in Fig. 2. The rescaled mass term is equal to 
zero on the critical line, thus the most relevant pertur- 
bation is the marginal energy-energy coupling [Eq. (25)]. 
Together with the action of the unperturbed Ising mo- 
dels from Eq. (9) the Ashkin-Teller field theory is readily 
identified:32'4i 



Aat = Aa + A% + dwdw eaeb 



(27) 



where /c oc is the renormalized coupling constant. The 
sign of the marginal perturbation can in principle matter, 
but for the present case it does not affect the conclusion 
that we stay on the critical line. Since the free fermionic 
theory is located in the middle of a line of c = 1 conformal 
field theories (see, for example, Ref. 26), the theory is 
critical on either side of the fixed point of two decoupled 
Ising models. 



D. Ashkin-Teller model 

The Ashkin-Teller lattice model was introduced as a 
generalization of the Ising model to a four-component 
system. 4^ However, in the 1970s it was shown that the 
model can be mapped onto a system of two Ising models 
{A and B) residing on the same lattice and interacting 
via an additional four-spin interaction:'*^ 



J4 E Si ^ • Sj ^ • Si ^ • Sj ^ 

NN 

(28) 



The field-theoretic action of this model is the same as in 
Eq. (27) with a coupling constant k — f{J,Ji). Thus, 
the four-spin interaction is mimicked by the perturba- 
tive onset of the energy-energy coupling, which is, in our 
case, proportional to the square of the original nearest- 
neighbor coupling Ji. The rich phase diagram of the 
Ashkin- Teller model given in coupling constants J and 



J4 in Refs. 2, 43, and 44 includes a critical line which 
represents the one-dimensional flow diagram of the cor- 
responding field theory [Eq. (27)] at criticality with a sin- 
gle parameter given by the renormalized coupling k. It 
starts for fc = at the conformally invariant fixed point 
of two decoupled Ising models and ends at the Potts- 
critical end point. In our calculations we arrived at a 
CFT with c = 1 in the presence of a marginal operator. 
Both are necessary conditions to fulfill so that the theory 
can exhibit varying critical exponents. ^'''^^ This scenario 
would be consistent with our numerical results for large 
J2 and earlier descriptions of the frustrated Ising model. ^ 
The Potts-critical end point would also allow for an onset 
of a non-critical line of first-order phase transitions as we 
see it in the MC simulations. However, since our effective 
field theory is only derived perturbatively, the critical be- 
havior does not necessarily need to be described by the 
Ashkin-Teller model up to Ji « J2. Thus, the explana- 
tion of the observed first-order transition (1 < J1/J2 < 2) 
might be beyond the scope of the Ashkin-Teller model. 



V. DISCUSSION 

Motivated by Refs. 13, 15-17 we had a closer look at 
the phase transition from the paramagnetic phase into 
the antiferromagnetic collinear ordered phase of the two- 
dimensional frustrated J1-J2 Ising model. We performed 
extensive additional MC simulations to verify the weak 
first-order character of the transition for the particular 
values of J2/J1 = 0.8 and 0.9 and recorded spin-spin cor- 
relation functions for larger values of J2 / Ji to gain an 
insight into the phase transition. However, the increas- 
ing length scales in the system do not allow for a reliable 
interpretation of the phase transition for Ji ^ J2 • Since 
for Ji = the model can be described by two decoupled 
Ising models and, hence, is treatable by means of CFT, 
we tried in a first attempt to calculate the corresponding 
central charge c with respect to the degree of frustration 
J1/J2 using transfer- matrix techniques. This computa- 
tion yielded the reasonable result c = 1 for small inter- 
coupling values Ji but suffered the same scaling problems 
as the MC simulations for larger values of J1/J2. Thus, 
we derived the continuous field theory for the discrete 
lattice model starting at the point Ji — and computed 
the perturbation induced by the nearest-neighbor cou- 
pling. Using OPE we closed the renormalization-group 
equations and arrived at an effective action which resem- 
bles the Ashkin-Teller model with c = 1. We have given 
symmetry arguments that highly relevant terms such as 
the spin-spin coupling operator can not be generated at 
any order. Thus, the critical behavior differs drastically 
from that of a two-layer Ising model. 

The derived effective field theory has a marginal 
energy-energy coupling as the most relevant perturba- 
tion to the decoupled Ising fixed-point action. Since 
this action equals the one of the Ashkin-Teller model we 
arrived at a c = 1 CFT which together with the pres- 
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ence of a marginal operator allows for continuous phase 
transitions with varying critical exponents. ^'^'^^ Thus, the 
weak-universality scenario of Landau and Binder^ is re- 
covered for Ji J2 but rather explained by the Ashkin- 
Teller model which opens the possibility for the occur- 
rence of the first-order phase transition wc observed in 
the MC simulations. The fact that the Ashkin- Teller 
model is only reached in the second order of the pertur- 
bation explains the large length scales that we find in 
the numerical analysis of the problem: the direct contin- 
uum limit yields only the irrelevant operators [Eq. (8)]. 
However, while these operators flow to zero under the 
renormalization group, they drive the marginal Ashkin- 
Teller coupling [Eq. (17)] to a finite value, thus generat- 
ing a crossover behavior as a function of the length scale. 
Moreover, it suggests that a further numerical analysis of 
the critical behavior would need very large systems sizes. 

It is also possible that the onset of a first-order transi- 
tion is caused by a level crossing of higher-energy states 
rather than by the renormalization fiow of the Ashkin- 
Teller model itself. 

Our results for the nature of the finite-temperature 



phase transition should also be relevant if small quantum 
fluctuations are included. Indeed Ref. 46 also finds a 
first-order transition close to the highly frustrated point 
which turns into a second-order transition for large J2 in 
a certain quantum generalization of the present model. 
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Appendix A: Second order perturbation using OPE 



The calculation of the second-order terms given in Eqn. (16)-(20) is shown excmplarily for the terms in Eq. (16) 
and (17): 



[Oi{z,z) - 02{z,z)] [Oi{w,w) - 02iw,w)] 
=Oi{z,z)Oi{w,w) ~ Oi{z,z)02iw,w) - 02{z,z)0i{w,w) + 02iz,z)02iw,iu) (Al) 
=dzaAiz,z)dzaBiz,z) dw<JA{w,w)dw<JB{w,w) - dzaAiz,z)d^aBiz,z) du,(JA{w,w)du,(TB{w,w) 

- dgaAiz,z)dgaBiz,z) dy,aA{w,w)dy,(jB{w,w) + dz(JA{zjz)dgaB{z,z) djBaAiw,w)du,<JBiw,w) . (A2) 



Using the fact that only products of operators on the same sublattice are non-zero and extracting the partial deriva- 
tives, we rewrite the last line before applying the general normalized OPE [given in Eq. (13)] for the spin-operator 
products: 



=dzdio<TA{z,z)aA{w,w) dzdw<7B{z, z)<7b{w,w) - dzdiaaA{z, z)aA{w,w) dzdwCrB{z, z)(jb{w,w) 

- dzdwCrA{z,z)aA{w,w) dzdiuaB{z,z)aB{w,'w)+dzdwaA{z,z)aA{w,w) dzdwCrB{z,z)aB{w,w) (A3) 

=dzd^, [\z - w\^^ + \\z ~ w\i£A{w,w) + \{z-w)^{z- w)i dn,£A{w,w) + H.c.^ 
X dzdw {\^z - w\^^ + \\z - w\^eB{w, w) + \{z - w)^ {z - ■w)idw£B{w, w) + H.c.^ 

- dzdw{- ■ ■A)dzdw{. . -b) - dzdw{. ■ ■A)dzdn,{. . .b) + dzdw{. . .A)dzdu,[. ■ -b) ■ (A4) 



In the following we only regard the regular terms and the pure energy terms for the sake of clarity. Note that 
\z — = {z — w)^ {z — w)^ which enters into the calculation of the partial derivatives. 

= ^ w)"T(z - + ^(z - (z - w)i£A{w, w) + O [3 ^^,^^£ a)) 

^ " - w)-i + -^{z ~ w)-^{z - w)ieBiw,w) + (A5) 

- 2 (-^(2 - w)-'i{z - + - wy^z ~ wyssAiw,-!!)) + 0{di^,a}eA)j 

^ ~ '^)"^(^ ~ '^)~^ + Y|g(z-w)~5(z-w)~t£B(u;,iy) + 0(a{^,c}e_B)^ (A6) 

+ ^ w)-i{z - w)-^ + -^{z - w)i{z - w)-^eAiw,w) + 0{d[^,^yeA)j 

{--ki^ -w)~Hz -w)-^ + ^{z -w)^z ~w)-^eBiw,w) + 0{d{^^.a}SB)'j ■ (A7) 
The calculation of the products (x) and summing over Eqs. (A5)-(A7) yields 

[Oi{z, z) - 02{z, z)] [Oi{w, w) - 02{W, w)] 
= ^ ((^ - w)~^{z - w)--^ + H.c.) -^{z~ - w)-^ 

- {-ik ((^ - w)-T(z - u,)i + H.c.) - ^{z - wyi{z - w)-i^ {ea + es) (A8) 

+ (liil +H.C.) - ^{z-w)-'i{z~wy^'^ EASE (A9) 

+ C(9{„,fi)}e7) 

which corresponds to the terms given in Eqs. (16) and (17). 
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